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Searches for dark matter imprints are one of the most active areas of current research. We focus 
here on light fields with mass ms, such as axions and axion-like candidates. Using perturbative 
techniques and full-blown nonlinear Numerical Relativity methods, we show that (i) dark matter 
can pile up in the center of stars, leading to configurations and geometries oscillating with frequency 
which is a multiple of / = 2.5 x 10^^ {rnB<? jeV') Hz. These configurations are stable throughout 
most of the parameter space, and arise out of credible mechanisms for dark-matter capture. Stars 
with bosonic cores may also develop in other theories with effective mass couplings, such as (massless) 
scalar-tensor theories. We also show that (ii) collapse of the host star to a black hole is avoided by 
efficient gravitational cooling mechanisms. 
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I. Introduction. The standard picture for the evolu¬ 
tion and structure of our universe relies on the existence 
of unseen forms of matter, generically called dark mat¬ 
ter (DM). The evidence for DM in observations is over¬ 
whelming, starting with galaxy rotation curves, gravi¬ 
tational leasing and the cosmic microwave background. 
While carefully concocted modified theories of gravity 
can perhaps explain almost all observations, the most at¬ 
tractive and accepted explanation lies in DM being com¬ 
posed mostly of cold, collisionless particles. 

Ultralight fields, such as axions or axion-like candi¬ 
dates are an attractive possibility [H-Q- Axions were 
originally devised to solve the strong-CP problem, but re¬ 
cently a plethora of other, even lighter fields with masses 
10-10 — 10-33 eV/c^, have also become an interesting pos¬ 
sibility in what is commonly known as the axiverse sce¬ 
nario [J. The simplest possible theory is that of a mas¬ 
sive scalar cj) or vector minimally coupled to gravity, 
and described by the Lagrangian 
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is the Maxwell 


tensor and F = The mass rriB of the boson un¬ 

der consideration is related to the mass parameter above 
through pLs,v = nT-B/fi, and the theory is controlled by 
the dimensionless coupling 
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where Mp is the total mass of the bosonic configuration. 

It is appropriate to emphasize that DM has not been 
seen nor detected through any of the known standard 
model interactions. The only evidence for DM is through 
its gravitational effect. Not surprisingly, the quest for 
DM is one of the most active fields of research of this cen¬ 
tury. Because DM interacts feebly with Standard Model 


particles, and thanks to the equivalence principle, the 
most promising channel to look for DM imprints consists 
of gravitational interactions. It turns out that classical 
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FIG. 1. Comparison between the total mass of a boson star 
(comp/es scalar) and an oscillaton {real scalar or vector fields), 
as a function of their radius R. R is defined as the radius 
containing 98% of the total mass. The procedure to find the 
diagram is outlined in the main text. 


massive bosonic fields minimally coupled to gravity, de¬ 
scribed by Eq. m - which will be our working model for 
DM here - can form structures Self-gravitating 

complex scalars may give rise to static, spherically- 
symmetric geometries called boson stars, while the field 
itself oscillates ii (for reviews, see Refs. [Tol - [l^ ). 
On the other hand, real scalars have a non-trivial time- 
dependent stress-tensor but may give rise to long-term 
stable oscillating geometries called oscillatons 0. Both 
solutions arise naturally as the end-state of gravitational 
collapse iHlil , and both structures share a similar 
mass-radius relation, summarized in Fig.[TJ In this figure 
we also show the mass-radius relation for massive vectors, 
which as far as we know has not been discussed before. 
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Boson stars and oscillatons have a maximum mass Mmax, 
given approximately by 


M„ 
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= 8 X 10"^^ 
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for scalars and slightly larger for vectors. Oscillatons 
are not truly periodic solutions of the field equations, 
as they decay through quantum and classical processes. 
However, their lifetime Tdecay is extremely large for the 
masses of interest mill, 
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II. Compact stars with DM cores. Although not 
framed in the context of DM capture, some stars with 
complex scalar fields in their interior were considered be¬ 
fore in the literature (Tsl - I^ . Our construction is more 
generic. 

Consider the Lagrangian (ED. augmented with the 
stress-energy tensor for a perfect fluid, 

TZd = ipF + P)u>^u’' + Pg>^\ (5) 


with the fluid four-velocity, 0,0, Oj, 

and pf, P its density and pressure. It is easy to 
show that the vector field must satisfy the constraint 
= 0, while the Bianchi identities and the 
Klein-Gordon equation, impose the conservation equa¬ 
tions = 0. Let us focus for simplicity on scalar 

fields. For complex scalars, self-gravitating solutions 
were considered many times, consisting of spherically 
symmetric, static boson stars [l0l - [T3j| . In Ref. Q it was 
shown that real massive scalars admit oscillating self- 
gravitating solutions. Consider therefore a general time- 
dependent, spherically symmetric geometry 

ds^ = -Fit, r)df + Bit, r)dr^ + r^dU'^ . (6) 

Rescaling the scalar field as (^ —>• and defining 

the function Cit,r) = Bit,r)/Fit,r), the field equations 
then lead to the system of partial differential equations 
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These equations suggest the following periodic expansion 

OO 
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where = (i?(t, r), C'(t, r), PFit, r). Pit, r)). The equa¬ 
tions of motion need to be supplemented by an equation 
of state. We will focus exclusively, and for simplicity, on 
a poly tropic equation of state P = Kp^ with AT = 100 
and r = 2, which can mimic neutron stars [20, but our 
results generalize to other equations of state Q. Inserting 
the expansion EID into the system ED- EHD and trun¬ 
cating the series at a given j 0, gives a set of ordinary 
differential equations for the radial Fourier components 
of the metric functions and the scalar field. These equa¬ 
tions need to be solved imposing regularity at the origin 
and asymptotic flatness. 

A useful quantity to describe scalar-fluid stars is the 
scalar field’s energy density, given by 

2p0 = —2 Tq® = Z'/ 9 tt + 4 ^ ! 9 rr + ■ (12) 

With this, we define the time-average total mass in the 
fluid and bosons as 

Mf,b= J (Zb pf,b^ r'^dr , (13) 

where <> denotes a temporal average. The total mass 
Mt can be found in the usual way through the met¬ 
ric component prr which asymptotically approaches the 
Schwarzschild solution at infinity. The procedure to find 
fluid-vector stars is similar and will not be described fur¬ 
ther. Figure |T] shows a phase diagram of solutions when 
there is no fluid. We find that the maximum mass for 
a stable scalar oscillaton is Mt ~ 0.6/ps in agreement 
with previous studies [3, [1^. As far as we know, the 
vector solution has not been described anywhere else. 

When the fluid is present, the solutions depend on two 
parameters pfo(O) and (/)i(0), and can be characterized 
by the mass coupling psMq, where Mq is the mass of 
the static star for vanishing scalar field, corresponding to 
the same value of central density pfo(O)- Different con¬ 
figurations are shown in Fig. |2l where we compare the 
scalar energy density against the fluid’s density. In all 
cases we also show the density profile when the scalar 
field is trivial, corresponding to the same fluid mass Mf- 
The overall behavior might have been anticipated from 
an analysis of Fig. |T] for light fields, psMq < 1, the 
scalar profile is extended and the pure oscillaton solution 
is broad and light. As such, the scalar has a negligible in¬ 
fluence on the fluid distribution (as can be seen from the 
fact that the zero-scalar line Ptj>=o overlaps with the fluid 
line), and these stars simply have an extended scalar con¬ 
densate protruding away from them. In fact, our results 
are compatible with a decoupling between the boson and 


^ In geometrical units G = c = 1 and F = 2, K has units length^ 
and can be used to set the length-scale of the problem in the 
absence of the scalar field [^ . The choice K = 100 was used in 
e.g. Ref. [531. 

^ We find that the solutions typically converge already for j = 2 
for most of the parameter space. 
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FIG. 2. Comparison between the (time average) energy 
density of the scalar field and the fluid pF for mixed 
scalar oscillatons and fermion fluids. From left to right, 
psMo = 0 . 1 , Mb/Mt ~ 21 %, psMq = 1 , Mb/Mt ~ 5%,and 
psMo — 10, Mb/Mt ~ 5%. Squares denote the correspond¬ 
ing quantities for complex fields (i.e. mixed boson-fluid stars). 
The overlap is nearly complete. Here Mo and P 0 =o are the 
total mass of the star and the energy density of the fluid, re¬ 
spectively, when the scalar field vanishes everywhere. In the 
left panel, the p</,=o and the pp lines are indistinguishable, 
because light fields have a negligible influence on the fluid 
distribution. 



FIG. 3. Density profile of the fluid at t = 0 and its first 
Fourier components for psMo = 10 and Mb/Mt ~ 5%. 


fluid for large pLs^B- For this case Fig.[T]alone is enough 
to interpret the bosonic distribution. For example, for 
PLsMq = 10, Mb/Mt = 5%, we get ~ 0.3, which 

would imply from Fig. [T] that psR ^ 20 for the scalar 
field distribution. This is indeed apparent from Fig. [2] 
Similar conclusions were reached when studying mixed 
fermion fluid/boson stars with complex fields E|. 
Global thermalization. We should stress that the com¬ 
posite stars we study here in general oscillate, driven 
by the scalar field. For these configurations, the con¬ 
vergence of the series m is rather fast, meaning that 


® In fact, the structure of mixed oscillatons and fluid stars is almost 
identical to that of geons and fluids, as can be seen from Fig. [2] 
where we overplot with dotted lines the complex field case. 


the amplitude of the oscillations is relatively small, but 
nonzero. For high scalar field central densities, first-order 
terms j = 1 in the density might become of the order of 
the zeroth-order term, making it difficult to find these 
configurations with good accuracy. For a given pFo{0), 
we expect these configurations to become unstable at 
some threshold c/i (0) > c/l (0). Our results are consistent 
with what was previously found for boson-fermion fluid 
stars (isl . [T^ . Although field configurations with high 
Mb/Mt are more challenging to find for large psMq, 
we expect them to follow the same kind of behavior as 
that found in boson-fermion stars. An example showing 
the amplitude of the oscillations is shown in Fig. for 
Mqps = 10, for which the star is oscillating with large 
amplitudes (the oscillating component is of the same or¬ 
der of magnitude as the static one). We find that even 
for Mb/Mt = 0.01, and for psMo = 10 the oscillations 
are of order 10% of the static component. These oscilla¬ 
tions imply that both the scalar field and the fluid density 
(which is coupled to it gravitationally), vary periodically 
with a frequency El 

/ = 2.5 X IQi^ Hz , (14) 

or multiples thereof. For axion-like particles with masses 
^ \0~^ eV/(?, these stars would emit in the microwave 
band. These oscillations are driven by the boson core 
and might have observable consequences; it is in principle 
even possible that resonances occur when the frequency 
of the scalar is equal to the oscillation frequency of the 
unperturbed star. The joint oscillation of the fluid and 
the boson might be called a global thermalization of the 
star, and is expected to occur also for boson-star-like 
cores (which give rise to static boson cores), once the 
scalar is allowed to have non-zero couplings with the star 
material. 

We have neglected viscosity in the star’s fluid and local 
thermalization. Viscous timescales for neutron star os¬ 
cillations can be shown to be large compared to the star 
dynamical timescale R, but small when compared to the 
(inverse of) the accretion rate likely to be found in any 
realistic configuration [2^. As such, we expect that vis¬ 
cosity will damp global oscillations of the star, eventually 
leading to a depletion of the scalar field core. A similar 
effect will occur with local thermalization of the scalar 
with the star material; detailed studies of these effects 
are still necessary. 

Stability. Finally, the stability of these solutions cannot 
easily be inferred from a mass versus radius diagram, and 
it requires a dynamical analysis which goes beyond the 
scope of this Letter. However, one can argue that a nec¬ 
essary condition for stability is that the binding energy 
Mt — Mp — Mb he negative M- Although negative val¬ 
ues do not necessarily imply stability, they do give strong 


The fundamental frequency a, ~ /rg v f^. Ilfill . 































4 


time: t/w=12.0 time: t/w=36.0 time: t/w=44.0 



FIG. 4. Snapshot of density p during the collision of two equal-mass oscillatons, each with Mps ^ 0.5, Rps 15, along the 
plane of collision, for half the space (the remaining space can be obtained by symmetry). Notice how the final configuration is 
more compact, and massive. At very late times, we find that it relaxes to a point on the curve of Fig. [T] 


support to the claim that these configurations are sta¬ 
ble. A more careful stability analysis shows that for suf¬ 
ficiently small </>! (0) and for stars which are stable in the 
absence of scalars, the negativityof the binding energy 
is a good criterion for stability (30|. To summarize, for 
sufficiently small scalars composite stars are dynamically 
stable. On the other hand, our results show that these 
configurations can be understood well from the mass- 
radius relation of oscillatons. The maximum mass sup¬ 
ported is ([3]), Mniax/Af© = 8 X 10~^^ eV/{mBC^), which 
for a neutron star and an axion field of mass 10“® eV 
falls well within the stability regime [3l|. 

III. Formation: DM accretion, core growth and 
gravitational cooling mechanism. We have shown 
that pulsating stars with DM cores exist as solutions of 
the field equations. Do they form dynamically Pulsat¬ 
ing purely bosonic states certainly do Bin, El- There 
are two different channels for formation of composite 
fluid/boson stars. One is through gravitational collapse 
in a bosonic environment, through which the star is born 
already with a DM core. The second process consists of 
capture and accretion of DM into the core of compact 
stars m, H^. 

Once the scalar is captured it will interact with the bo¬ 
son core. Interactions between complex fields have shown 
that equal mass collisions at low energies form a bound 
configuration [s^, H^. In other words, two bosonic cores 
composed of complex fields interact and form a more 
massive core at the center. This new bound configu¬ 
ration is in general asymmetric and will decay on large 
timescales [T^, the final state being spherical [s^. We 
have repeated this analysis for two oscillatons colliding 
at sufficiently small energies, using the methods of Nu¬ 
merical Relativity outlined in Refs. [H, [s^. We find 
that generically the collision of two oscillatons whose to¬ 
tal mass is below the peak value will give rise to a more 
compact, more massive oscillaton configuration. By con¬ 
trast, the collision of oscillatons whose mass exceeds the 
peak value results in large mass losses and a lighter final 
configuration. An example of such a collision is shown in 
Fig. Ufor equal-mass oscillatons each with Mpg-s ~ 0.5 
(close to the peak value of Fig. [I}. Note that the total 
mass is larger than the peak value and one would naively 


predict gravitational collapse to a black hole (BH). In¬ 
stead, the collision first produces a massive, more com¬ 
pact oscillating object, accompanied by large scalar field 
losses, as is apparent from Fig. |4l On large timescales 
the configuration relaxes to a very low-mass structure. 

This is a general feature of what has been termed 
the “gravitational cooling mechanism:” a very efficient 
(dissipationless) mechanism that stops them from grow¬ 
ing past the unstable point, through the ejection of 
mass [l^, . Such features have been observed in the 

past in other setups, such as spheric ally symmetric grav¬ 
itational collapse (see Fig. 2 in Ref. [3), or slightly per¬ 
turbed oscillatons [l^. Gravitational cooling provides a 
counter-example to an often-used assumption in the liter¬ 
ature: that stars accreting DM will grow past the Chan¬ 
drasekhar limit for the DM core and will collapse to a 
BH [H, . Our results show that this need not be 

the case, if the DM core is prevented from growing by a 
self-regulatory mechanism, such as gravitational cooling. 

Thus, even though other more detailed simulations are 
still needed, the likely scenario for evolution would com¬ 
prise a core growth through minor mergers, slowing down 
close to the mass-radius peak (see Fig. ^ at which point 
it stops absorbing any extra bosons 1^, [2^. In other 
words, the unstable branch is never reached. This phe¬ 
nomenology is specially interesting, as it would also pro¬ 
vide a capture mechanism for these fields which is inde¬ 
pendent of any putative nucleon-axion interaction cross 
section: as we discussed, the bosonic core grows (in mass) 
through accretion until its peak value. At its maximum, 
it has a size Rb/Mq ^ Mb/Mq. This is the bosonic 
core minimum size, as described by Fig. [TJ Even for 
Mb = O.OIMt the core has a non-negligible size and is 
able to capture and trap other low-energy oscillatons. 
IV. Conclusions. Previous works on the subject of DM 
accretion by stars have implicitly assumed that the DM 
core is able to grow without bound and eventually col¬ 
lapse to a BH Our results, from full nonlinear 

simulations of the field equations, show that the core may 
stop growing when it reaches a peak value, at the thresh¬ 
old of stability, if DM is composed of light massive fields. 
Gravitational cooling quenches the core growth for mas¬ 
sive cores and the core growth halts, close to the peak 
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value (c.f. Fig. [T|). 

Finally, our results are formally valid only for the the¬ 
ory o, and require the fields to be massive. Minimally 
coupled, multiple (real) scalars, interacting only gravi¬ 
tationally, were also shown to give rise to similar con¬ 
figurations [ 4 ^. We have explicitly verified that similar 
results may hold in other setups. For example, scalar- 
tensor theories (which can be shown t o g ive rise to an 
effective position-dependent mass term |44l.l45l|l and non- 
minimally coupled massless scalars are, in principle, also 
able to develop nontrivial configurations with a scalarized 
core. Further details will be presented in a forthcoming 
work. 
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